FITZHUGH_NAGUMO

Overview

The FITZHUGH_NAGUMO function numerically solves the FitzHugh-Nagumo (FHN) system of ordinary differential equations, a simplified model of excitable systems such as neurons. The FitzHugh-Nagumo model was developed by Richard FitzHugh in 1961 and Jinichi Nagumo et al. in 1962 as a two-dimensional simplification of the more complex Hodgkin-Huxley model.

The FHN model is a relaxation oscillator that captures the essential dynamics of neuronal spike generation. When an external stimulus exceeds a threshold, the system exhibits a characteristic excursion in phase space before returning to rest. The variable V represents the membrane potential (fast activation), while W is a recovery variable (slow inactivation) representing the combined effect of sodium channel reactivation and potassium channel deactivation.

The system is governed by the following coupled differential equations:

\frac{dV}{dt} = V - \frac{V^3}{3} - W + I_{\text{ext}}
\frac{dW}{dt} = \frac{1}{\tau}(V + a - bW)

where a and b are model parameters controlling the recovery dynamics, \tau is a time-scale separation parameter, and I_{\text{ext}} is the external applied current.

This implementation uses SciPy’s solve_ivp function to perform numerical integration. The solver supports multiple integration methods including RK45 (default), RK23, DOP853, Radau, BDF, and LSODA. Explicit Runge-Kutta methods (RK45, RK23, DOP853) are generally suitable for non-stiff problems, while implicit methods (Radau, BDF) handle stiff systems more efficiently.

The function returns time series data for both variables, enabling analysis of action potential dynamics, limit cycles, and bifurcation behavior under varying stimulus conditions.

This example function is provided as-is without any representation of accuracy.

Excel Usage

=FITZHUGH_NAGUMO(v_initial, w_initial, a, b, tau, i_ext, t_start, t_end, timesteps, solve_ivp_method)
  • v_initial (float, required): Initial membrane potential.
  • w_initial (float, required): Initial recovery variable.
  • a (float, required): Model parameter a in the recovery equation.
  • b (float, required): Model parameter b in the recovery equation.
  • tau (float, required): Time scale parameter for the recovery variable.
  • i_ext (float, required): External current applied to the neuron.
  • t_start (float, required): Start time for integration.
  • t_end (float, required): End time for integration.
  • timesteps (int, optional, default: 10): Number of timesteps to solve for.
  • solve_ivp_method (str, optional, default: “RK45”): Integration method for solve_ivp.

Returns (list[list]): 2D list with header [t, V, W], or error message string.

Examples

Example 1: Demo case 1

Inputs:

v_initial w_initial a b tau i_ext t_start t_end
0 0 0.7 0.8 12.5 0.5 0 1

Excel formula:

=FITZHUGH_NAGUMO(0, 0, 0.7, 0.8, 12.5, 0.5, 0, 1)

Expected output:

t V W
0 0 0
0.1111 0.05839 0.006455
0.2222 0.1228 0.01341
0.3333 0.1938 0.02091
0.4444 0.272 0.02902
0.5556 0.3576 0.0378
0.6667 0.4509 0.04731
0.7778 0.5515 0.0576
0.8889 0.659 0.06875
1 0.772 0.0808

Example 2: Demo case 2

Inputs:

v_initial w_initial a b tau i_ext t_start t_end timesteps solve_ivp_method
0 0 0.7 0.8 12.5 0.5 0 1 10 RK23

Excel formula:

=FITZHUGH_NAGUMO(0, 0, 0.7, 0.8, 12.5, 0.5, 0, 1, 10, "RK23")

Expected output:

t V W
0 0 0
0.1111 0.05839 0.006455
0.2222 0.1228 0.01341
0.3333 0.1938 0.02091
0.4444 0.2719 0.02901
0.5556 0.3575 0.03779
0.6667 0.4507 0.04729
0.7778 0.5513 0.05758
0.8889 0.6588 0.06873
1 0.7718 0.08078

Example 3: Demo case 3

Inputs:

v_initial w_initial a b tau i_ext t_start t_end timesteps solve_ivp_method
1 0.5 0.7 0.8 12.5 0.5 0 2 10 RK45

Excel formula:

=FITZHUGH_NAGUMO(1, 0.5, 0.7, 0.8, 12.5, 0.5, 0, 2, 10, "RK45")

Expected output:

t V W
0 1 0.5
0.2222 1.144 0.5242
0.4444 1.272 0.5506
0.6667 1.379 0.5786
0.8889 1.462 0.6079
1.111 1.52 0.638
1.333 1.558 0.6685
1.556 1.579 0.6991
1.778 1.589 0.7296
2 1.59 0.7597

Example 4: Demo case 4

Inputs:

v_initial w_initial a b tau i_ext t_start t_end timesteps
-1 0.2 0.5 0.7 10 0.3 0 1.5 10

Excel formula:

=FITZHUGH_NAGUMO(-1, 0.2, 0.5, 0.7, 10, 0.3, 0, 1.5, 10)

Expected output:

t V W
0 -1 0.2
0.1667 -1.093 0.1886
0.3333 -1.181 0.1758
0.5 -1.261 0.1618
0.6667 -1.332 0.1467
0.8333 -1.391 0.1308
1 -1.44 0.1141
1.167 -1.479 0.09685
1.333 -1.507 0.07927
1.5 -1.527 0.06148

Python Code

from scipy import integrate as scipy_integrate
import numpy as np

def fitzhugh_nagumo(v_initial, w_initial, a, b, tau, i_ext, t_start, t_end, timesteps=10, solve_ivp_method='RK45'):
    """
    Numerically solves the FitzHugh-Nagumo system of ordinary differential equations for neuron action potentials using scipy.integrate.solve_ivp.

    See: https://docs.scipy.org/doc/scipy/reference/generated/scipy.integrate.solve_ivp.html

    This example function is provided as-is without any representation of accuracy.

    Args:
        v_initial (float): Initial membrane potential.
        w_initial (float): Initial recovery variable.
        a (float): Model parameter a in the recovery equation.
        b (float): Model parameter b in the recovery equation.
        tau (float): Time scale parameter for the recovery variable.
        i_ext (float): External current applied to the neuron.
        t_start (float): Start time for integration.
        t_end (float): End time for integration.
        timesteps (int, optional): Number of timesteps to solve for. Default is 10.
        solve_ivp_method (str, optional): Integration method for solve_ivp. Valid options: RK45, RK23, DOP853, Radau, BDF, LSODA. Default is 'RK45'.

    Returns:
        list[list]: 2D list with header [t, V, W], or error message string.
    """
    # Validate input types
    try:
        v0 = float(v_initial)
        w0 = float(w_initial)
        a_param = float(a)
        b_param = float(b)
        tau_param = float(tau)
        i_ext_param = float(i_ext)
        t0 = float(t_start)
        t1 = float(t_end)
        ntp = int(timesteps)
    except Exception:
        return "Invalid input: All initial values, parameters, and timesteps must be numbers."
    if t1 <= t0:
        return "Invalid input: t_end must be greater than t_start."
    if tau_param <= 0:
        return "Invalid input: tau must be positive."
    if ntp <= 0:
        return "Invalid input: timesteps must be a positive integer."
    valid_methods = ['RK45', 'RK23', 'DOP853', 'Radau', 'BDF', 'LSODA']
    if solve_ivp_method not in valid_methods:
        return f"Invalid input: solve_ivp_method must be one of {valid_methods}."
    # Create time array for evaluation
    t_eval = np.linspace(t0, t1, ntp)

    # FitzHugh-Nagumo ODE system
    def fhn_ode(t, y):
        v, w = y
        dv_dt = v - v**3 / 3 - w + i_ext_param
        dw_dt = (v + a_param - b_param * w) / tau_param
        return [dv_dt, dw_dt]
    # Integrate
    try:
        sol = scipy_integrate.solve_ivp(
            fhn_ode,
            [t0, t1],
            [v0, w0],
            method=solve_ivp_method,
            dense_output=False,
            t_eval=t_eval
        )
    except Exception as e:
        return f"Integration error: {e}"
    if not sol.success:
        return f"Integration failed: {sol.message}"

    # Helper to check for NaN or Inf
    def is_finite(val):
        return val == val and abs(val) != float('inf')

    # Format output: header row then data rows
    result = [["t", "V", "W"]]
    for i in range(len(sol.t)):
        t_val = float(sol.t[i])
        v_val = float(sol.y[0][i])
        w_val = float(sol.y[1][i])
        # Disallow NaN/Inf
        if not (is_finite(t_val) and is_finite(v_val) and is_finite(w_val)):
            return "Invalid output: NaN or Inf detected in solution."
        result.append([t_val, v_val, w_val])
    return result

Online Calculator